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We present a microscopic approach to the quantum tunnehng of vortices. The formalism charac- 
terizes the rate at which a many-body superconducting state with a vortex in one location makes 
a transition to a second many-body superconducting state with a vortex in a second location. The 
theory is utilized to study the effect of a d-wave order parameter upon the directionality of vortex 
tunneling. 
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For decades, researchers have studied the motion of 
vortices in type- II superconductors [l|. This motion de- 
termines fundamental characteristics of superconductors, 
including critical current and mixed-state resistivity [3, 
which often dictate their technological usefulness. Be- 
cause of the complicated nature of a vortex in a super- 
conductor, however, many aspects of the theoretical un- 
derstanding of vortex motion have remained incomplete. 
In particular, the theory of quantum tunneling of vortices 
|3,|j| has been phenomenological in character, treating a 
vortex as a point-particle governed by an effective ac- 
tion. In this paper, we approach the problem of vor- 
tex quantum tunneling from a fully microscopic quan- 
tum mechanical perspective. Our calculations are based 
on the Nambu-Gorkov equations of superconductivity; 
we do not require effective parameters and do not as- 
sume that the vortex can be treated as a point-particle. 
Because of its microscopic foundation, the theory (i) is 
capable of addressing qualitative phenomena that are be- 
yond the reach of effective models and (ii) provides a 
foundation for quantitative, first principles analysis of 
vortex motion in real materials. We focus in this pa- 
per upon non-dissipative motion, a topic that has elicited 
sustained theoretical concern even without the additional 
complication of dissipation 5, 6] (dissipation would be in- 
corporated into the theory in the spirit of 7] ) . As a first 
application, we study the effects of the d-wave nature of 
the order parameter in high-temperature superconduc- 
tors |S|. Our calculations predict that d-wave symmetry 
suppresses tunneling along the nodes of the order param- 
eter. This effect cannot be captured even in principle 
within the existing point-particle theories |9|- 

We begin our formulation by introducing the state \vn) 
that describes a superconductor with N electrons form- 
ing a vortex at the origin. For simplicity, we assume that 
the electrons are confined to 2 dimensions, which is ap- 
propriate for a superconducting thin film or possibly for 
a layered superconductor. The extension to 3 dimensions 
is straightforward. The microscopic Hamiltonian of the 
system is H — Ho + V where Ho describes electrons in an 
ionic potential and a homogeneous magnetic field, under- 



going an effective electron-electron interaction that gives 
rise to superconductivity. A separate pinning potential 
V results from microscopic impurities or dislocations. 

We consider the dynamics of this Hamiltonian in a ba- 
sis of states of the form \vn{p)) = T{p) \vn), where the 
operator T{p) translates the vortex from the origin to the 
position p. It is important to emphasize that we are not 
making a phenomenological point-particle |2| assump- 
tion about the vortex coordinate p when we deal with 
states of the form |wjv(p)); the vortex states are the low- 
est energy states of a superconductor in a magnetic field 
|10l | and simply compose the low-energy Hilbert space 
of the microscopic Hamiltonian. Physically, this basis 
is a natural choice for describing a superconductor with 
a moving vortex. To compute the behavior of a time- 
dependent vortex state \^{t)) = ^p,'^{p' ,t)\vN{p'))j 
in the usual way we first compute the eigen- 
states IV') — X]p' V'(p') kiv(p')) of the time-independent 
Schrodinger equation Ho + V^ — i? |^) = 0. In our basis, 
this equation is 

^ {vn{p)\Ho + V-E \vn{p')) i>{p') = . (1) 
P' 

Before we can contemplate solving this equation, the 
essential challenge is to evaluate the matrix elements 
{vn{p)\Ho + V\vn{p')) and the overlap {vn{p) \v]si(,p')), 
keeping in mind that the states in our basis are not or- 
thogonal. These matrix elements seem to require 2N 
dimensional many-body integrals, which can become ut- 
terly intractable for even relatively small TV. In the fol- 
lowing, we provide tractable forms for them by describing 
\^n{p)) using the mean-field, Nambu-Gorkov equations 
of superconductivity [ll|. Of course, a solution to (^ 
is a superposition 1-0) = X]p' ^{p') \^n{p')) of mean-field 
states, so that our analysis goes well beyond mean-field 
theory. 

To employ Nambu-Gorkov theory, we introduce the 
state \v) of a superconductor that has a vortex at the 
origin and is in contact with a particle bath. This 
state can be well characterized using mean-field theory 
and the N electron state \vn) can be obtained from 



it by projection \vn} = Pn\v)- The mean-field de- 
scription of \v) is contained in the Green's functions, 
G,^s'{r,r',t) = -{v\T{d,^,it),dl,^^,{t' = 0)} \v) which 
are defined using a particle-hole transformation from the 
usual electron destruction operators to dr,s=i = Cr,s=i 

After Fourier transforming in 



and (ir,s=-i 
time, one 



finds ^ 
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where K — 7r—(r) 



'^A)^ + U{r) - Ef and ^ is a 2 by 
2 matrix with elements G's.s'(r,r', i?). The crystal po- 
tential f/(r) includes the ionic potential and electron- 
electron Hartree interactions. The superconducting order 
parameter operator A is defined by A Gs.s'{i,J^' ,E) = 
J di: A(r,f) Gs,s'{i,J^' ,E) and should be computed self- 
consistently. This equation can be solved numerically 
on a finite real space lattice, in which case the Green's 
function can be regarded as a matrix. 

The vortex state |wAr(p)) is obtained by translating 
\vn), which contains a vortex at the origin. Because 
of the presence of an external, homogeneous magnetic 
field in the Hamiltonian, one effects the translations 
using magnetic translation operators [12, I13. Mag- 
netic translation operators are usually defined with re- 
spect to single-particle wavefunctions; in the Landau 
gauge A(r) = —Byx, we have the definition r(p)-0(r) = 
^z{7rB/^„)p^ix~pj2)^(^^ _ p^ ^j^gj.g ^^ ^ hc/2\e\ is the 

superconducting flux quantum I13. This operator must 
be generalized in order to translate a many-body vortex 
state properly. A suitable generalization is 



In this expression we have introduced fermionic operators 
of the form C^ that create electrons in an artificial, aux- 
iliary space. Such operators are needed because "hold- 
ing" space is required when exchanging a collection of ob- 
jects - to exchange the balls in two boxes, a third "hold- 
ing" box is required temporarily. In other words, one 
should not shift an electron from r to r -|- p without first 
placing an electron already in r -I- p into a holding box. 
One can verify that T(p)4^_p ,^4^_p ,^...c|;^_p^^_^ = 

gi(TrB/*<,)p„(a;i-p^/2)^t _^gi{-^B /<S>a)py(x„~p^ /2) ^f ^^g 

expected. If needed, the number projection operator Pjv 
can be conveniently incorporated into T[p) as follows 






After applying the particle-hole transformation to this 
operator, one can use Wick's theorem |lj| to eval- 
uate the overlap matrix element {vn{p) l^Nip')) — 
{v\T\p)T{p')Pm\v) that appears in equation ^. Up 
to an unimportant phase factor that can be absorbed 
into the vortex state, the result is 



{vn{p) \vn{p')) - / ge'(^°/2-^)xdet M{x). (3) 



KP) R.sl R,S R^s c-^^s ^,S) where iVo is the number of single-particle orbitals in the 

X nr.,(cr-p,.ct_p^^ + e'(-^/*'')^«(--''^/2)ct ,cr-p^,). systcm and 
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r,s);(r',s') 



-Gs^s'{r - p', r' -p',t = 0+) + G,,,,(r - p', v' p,t = Q-y''{x+{^B/'i>^)(py'(x'-pJ /2)-py(x' -p^/2))) 

I 



is an No x No matrix indexed by row (r, s) and column 
(r',.'). 

Note that equation ^ involves only determinants and 
a single dimensional integral. In place of a totally 
unmanagable 2N dimensional many-body integral for 
{vn[p) \vn(,p')), it provides a computationally tractable 
expression that we evaluate below for N up to 312. A 
clear intuitive interpretation of lO can be obtained by 
analogy to a simpler case. Imagine that \vn) were, in- 
stead of a superconducting vortex state, a Hartree-Fock 
state of N spinless fermions in a Hilbert space of total 
dimension N^ . In that case, \vm) would be a Slater 



determinant of N orbitals, C,i{v) (i — 1,...,N) drawn 
from a complete, orthonormal basis of No orbitals, Q (r) 
{i — l,...,No)- It would follow from the definition of 
the Slater determinant that {vn{p) \vn{p')) = det Z , 
for Z an A^ X iV matrix of overlap integrals Zi,j = 
J dr(* {r — p)(j {r — p') . Defining a related NoXNo matrix 



5,.. 



-'^,3 



p)Uv-p') ] = l,...,N 

j^N + l,...,No 



we find that detZ 
functions G'(r, r',t 



G(r,r',t = 0^ 



= detZ. Using the Green's 

0-) = E.<jvC.(r)C(r') and 

Si>Af C»(r)C(r'): we construct the 



^N, 



Mr,r' = -G(r-p',r'-p',t = 0+)+G{r-p' ,r'-p,t = 0") 

This matrix satisfies U'^ MIA — Z for the unitary matrix 
Ur,i = Cj(i" - p') ■ Therefore, {vn{p) \vn{p')) = dot Z = 
det Z — det A4. The equation is reasonable: the G(r — 
p', r' — p,t = 0~) part of Ai translates the N occupied 
orbitals in |ujv(p')) to p while the —G{r — p' , r' — p' ,t = 
0+) = d{r - r') - G{r - p' , r' - p',t ^ 0") part leaves 
the No — N unoccupied orbitals untranslated. Similar 
intuition applies to our superconducting case, although 
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we have a magnetic field and we also integrate over phase 
X to project \v) on to the state \vn). 

The reaso ning leading to equation Q also permits us 
to compute [lj| the matrix elements of the pinning po- 
tential 



{vn{p)\V\vn{p')) - {v\tHp)VT{p')Pn\v) (4) 

27r 



= E V.Aro) J ^e^'''"r^-''>''dctC{x,ro,So) 



where 



'C(x,ro, So)(r,s);(r',s') = -Gs.s'{r - p',r' - p',t = 0+)(l - S(^r',s');{r„,s„)SsaS) 

+Gs^s'{r - p', r' - p, i = o-)e^^'i^+i-^/'^^'>iPy'(-'-P^'/^)-Py(-'-P^/^y^) (l - <5(r^v);(r„,s„)'5.„,-l) . 

I 



Equations (P), ||3J), and 10} provide the ingredients for 
a microscopic calculation of quantum vortex tunneling. 
We stress that, although formulae (O and ^ are compli- 
cated in appearance, they are attractive and manageable 
from a computational standpoint. No formula is needed 
for the remaining matrix element (ujv(p)| Ho \vn{p')) if 
we are willing to approximate that Ho has |uAr(p)) as 
an eigenstate. We then just get {vm{p)\Ho\vn{p')) = 
Ev (vNip) \vn{p')) and the constant energy Ey merely 
shifts the eigenvalue i? in (^). 

This formalism makes is possible to study, for in- 
stance, how d^2_y2-wa,ye order parameter symmetry in- 
fluences vortex quantum tunneling in a high temperature 
superconductor. We model a high temperature super- 
conductor with the Hamiltonian described in reference 
|l5| . We follow the methodology of [l5j to obtain self- 
consistent quasiparticle amplitudes (u, v) and then sum 
these amplitudes to produce the self-consistent Nambu- 
Gorkov Green's function appearing in J^J) and Q) - this 
is equivalent to solving the Nambu-Gorkov equations di- 
rectly. Calculations are made on a 14 x 28 real-space 
lattice {No = 2 x 14 x 28 = 784 single particle states), 
where the number N « 0.8(A^o/2) = 312 of electrons 
in the superconductor is set near optimal doping. Fig- 
ure (Q) plots the magnitude of the calculated overlap 
1 {vn{p = Xx -\- Yy) \vm{'^)) \ as a function of displace- 
ment X and Y in units of the lattice constant a. Results 
are shown for an d2,2_j^2-wave superconductor; because it 
is a logarithmic plot, the figure is not noticeably changed 
for s-wave order parameter calculations. The results pre- 
sented in figure (^ do not depend sensitively upon the 
14 X 28 lattice size used in the numerical calculation. To 
demonstrate this, figure ^ displays the nearest neigh- 
bor overlap | {vn{p = ay) |wjv(0} | for lattices of increas- 
ing size. Note that the magnitude of the overlap changes 
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FIG. 1: Magnitude of the overlap 

\{vm{p = Xx + Yy) |ujv(0)) I of two vortex states as a 
function of displacement Xx + Yy between the states. 



little even when the number of electrons grows by more 
than a factor of 3. This finding confirms that interac- 
tions can lead to finite vortex overlaps in infinite systems 
[3 . The result is in contrast to the case of a system of N 
non-interacting bosons, for instance, in which the overlap 
decays exponentially with increasing number of particles. 
We model quantum vortex tunneling between two pin- 
ning sites by inserting impurities into our self-consistent 
calculation of the Green's function. Of course, if our aim 
were to guide technological efforts to raise critical cur- 
rents, we would seek quantitative results from the for- 
malism by computing a realistic pinning potential V . 
To explore generic behavior, though, we take the pin- 
ning potential to be a double well of well spacing p^, 
with depth given by 1/7 the Fermi energy, a reason- 
able magnitude. We further assume that the effect of 
this potential on the system's electrons is to lower the 
energy of \vn{0)) and \vn{Po)) below the energy of all 
other states \vn{p))- Then the Schrodinger equation (Q 
is solved, in the two dimensional Hilbert space spanned 
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FIG. 2: Magnitude of nearest neighbor overlap 

I {vn{p ~ ay) ji'jv(O)) I for calculations on lattices of in- 
creasing size, with s-wave and d^2_y2 wave order parameters. 



by these states. If the superconductor begins in a time- 
dependent state |*(i)) satisfying |\l/(t = 0)) = |ujv(Po))' 
the probability that tunneling has occured by time t is 
P(t) = 1 - \{vNiPo) l*(i))P • This probability is found 
to be periodic in our double well model, and we there- 
fore set the tunneling rate T{p^) equal to the reciprocal of 
the period. Although the model produces tunneling back 
to the initial well, it is assumed, as in other macroscopic 
quantum tunneling models [l£j , that the environment in- 
tercedes to prevent this from occuring. Since the dissipa- 
tive effects of the environment are not explicity included 
in the calculation, we should obtain an upper bound on 
the tunneling rate [17|. One could include the effects of 
dissipation in future work by recasting the dynamics of 
(E3) in a functional integral as in Q . 

As a measure of the directionality of the tunneling, 
we compute the ratio T{ax + ay)/T(ay), using formulae 
(P), 0, and 10} to evaluate T{pg). In the s-wave case 
for a 14 X 28 lattice, we obtain a ratio of 50%. In the 
dx^-y2-wave case for a 14 x 28 lattice, the ratio is 2%. 
In comparison to the s-wave case, the dx^^y2--wa,ve or- 
der parameter symmetry strongly suppresses the relative 
tunneling rate along the order parameter nodes. When 
the high potential region between the two wells is ob- 
structing motion along the nodes of the order parameter, 
it is more effective in confining the superconducting state. 

The effect could be probed by magnetization measure- 
ments of the kind reviewed in [l9j. A potential exper- 
imental arrangement involves the fabrication of two su- 
perconducting cylindrical tubes of square cross-section. 
One cylinder has the nodes of the d-wave order parame- 
ter parallel to its walls; the other cylinder has the nodes 
angled at 45° to its walls. A vortex trapped in the hollow 
interior of a cylinder by field cooling will tunnel out at a 
rate influenced by the order parameter orientation. 

The author thanks A. Auerbach, who initially sug- 
gested that the author study anisotropic dynamics of d- 
wave vortices, J. Jain, G. Koren, M. Wolraich and Y. 
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Space constraints make it impossible to present detailed derivations of our results within the body of our manuscript. 
This appendix supplies explicit mathematical arguments for the convenience of interested readers. 

PARTICLE HOLE TRANSFORMATION 

cn", 

^^ I To evaluate the overlap matrix element (wAr(p) \vn{p')), one inserts the appropriate translation operators 

^^ ■ {vNip) \vNip')) — {v\T^{p)T{p')Pn\v). Note that there is no operator P^r next to T\p); a second projection 

operator would not change the matrix element. 
C^H' Now, the definition of the many-body translation operator is 

(N 



X 



and it is convenient introduce a particle-hole transformation. We write 



o 
o 

Q^- = nRCR,_inR(c^_iCR,_i + CR^_iC^_j 

d . The phase factor at the end arises since magnetic boundary conditions imply nrCr-p.-ie*'^'^'^/*"^''!'*^^^''^"/^-' = 

^ [ nrCr,-ie*^'^-^/^*°^''''^"=. Defining "holding-space" operators -Dr,i = Cr,i and P'r,-! = Cj^ _^, we therefore find 

'^ : Tip) = nre*(-^/2*°)''^''''(nRi?J,^_jnR,5pR.s<5 + <s^R,5) 

The leading phase factor can be absorbed into the definition of the vortex state if desired, as \vn{p)) = 

(^ ' Y[^e~^^^^/'^°'>P''Py/'^T{p) \vn), and can therefore be omitted from (QJ. 

J> ■ A nearly identical calculation yields the operator T{p)Ppf, although an additional phase factor e*'-^"/^-'^ arises. 

^^ ' Knowing the translation operator, we can compute the overlap and translation matrix elements. 

m ■ 

o 



OVERLAP MATRIX ELEMENT 



The overlap matrix element is 

"c^ : {vn{p) \vn{p')) 

^ : ^ii,s{Dn^sDl,s + ^R,srfR,s)(nRi^R,-i) X iU^'Dl,_^)UR,^s'{Dn,^s'Dl,s, + 4',S'^R',s') 

.. i 1-5 \ r'—p ,s' r' — p,s' ' r'—p ,s' ■' ^ » » , j-/ i / \ / 

27r L 

-Gs,s'ir^p',r'-p',t^O+)] 

%^(^^/^-^)xdctMix). (3) 

ZTT 

The equality of (01 and ^ is demonstrated as follows. We start by noticing that the middle factors 
(nRZ?R._i)(nR'Z3j{^, _i) give unity when acting to the left or the right. With these operators eliminated, we consider 
the new middle factors 

Ur,s{Dr,sDIs + ^LsC^R,s)^R^s'(^R^5'^k^s' + 4',s'^R',s') = 

nR,s(£'R,s£'Ls + c^Ls^R^s + dIiA^^s + Dlf^,sDR,sdR,sdl,^s) -> 1 



since the operator d^s yields zero when it acts to the right, the operator d^ g yields zero when it acts to the left, 
and the combination dR^sc^R 5 yields unity. This leaves 



P 

Expanding the product out as a sum, one finds 

N, 



2tt 

27r ' ' H' 



I '^e,KNo/2-N)xy^ y^ ^^1 /'g-J5l(7rB/*„)p„(xi-p,/2)^t _ 

n— ri,si:...;rn,Sn 

rn-p,SnJ y ri-p ,si ' ■ ' 

rri-p ,s„J y r„+i-p ,s„ + i r„ + i-p ,s„ + i rN„-p,sjv„ tn^-P ,SNaJ 



Invoking Wick's theorem, this leaves us with a sum over permutations P 



/"^gi(W„/2-iV)xy^('_;^^Py^ y^ fg»^l(x+(^B/*o)(p;(a:i-p;/2)-p„(a;i-p,/2))) , .^t ^ , |^\ 

P n— ri,si;...;rTi,s^ 

g*s„(x+(7rB/$o)(p;,(a;„-p^/2)-p„(a:„-p,/2))) / i^t ^ , |^\A 



{v\dr „,.spd'' , \v)...{v\d _ / d^ , \v) 

= /" ^„i(W„/2-Ar) 



g.(W„/2-iV)x^(_l)Pn, (e-.(x+(-B/*„)(p;(.„-pL/2)-p„(.„-p./2))) {v\dl_p^^4^^ _p, ^^ \v) 



27r 

\ I rc.-rt ,sp. n-p ,Si ' ' y 



rp^ -P 



27r L 

-a,.'(r-p',r'-p',t = 0+)] 

^g.(w„/2-iv)xdetM(x) 

ZTT 

The symbol (— f)'^ takes the value +1 for even permutations and —1 for odd permutations. In the second to last 
equality, we recognize that a sum over permutations can be written as a determinant. 

MATRIX ELEMENT OF THE POTENTIAL 

The analysis that enable us to compute {vn(p) \vn{p')) can be extended to obtain matrix elements of the pinning 
potential (w^(p)| V \vNip')) = {v\ T\p)VT{p')Pn \v). 

First, the potential is expressed in terms of the d operators 



Given the form Q, it is necessary to compute {vn{p)\ dl ^ dr^.s„ \vn{p'))- The result is 
{vn{p)\ dl^^, Jr^^s, \vn{p')) 
= J ge^(^°/^-^)'^ {v\ {uM,^,)IlrAdr-pA^p.s + e--(-^/*°)''«(--''^/^)4-p.^r,.) 

= / ge^(^°/2-^)'^ {V\ (Hri^t _l)nr,.(l + e-(-^/*°)''»(^-''-/^)4-p,.^r..)i?i„^.„i?r„,.„ 
Ilr',Adr'-p'M'^p',s' + e'^'("+^^''/*°'''"^^'-''^/'')i^^,.-dr'-p',.')(nr'^r',-l) b) 

From this calculation, we find that 

{vn{p)\V\vn{p')) 

To, So P 

-L\ t^/ i-i,ai I \ I r^ p,Sj Yp^—p ,sp^ ' ' ^ ' rp. — p ,sp^ Ti—p.Si I 

\ To— p, — i vp^—p ,sp^ ' 'y 

+Fi(r„) (^e<'^+('^^/*°)K(-o-P;/2)-P«(-o-P./2))) (^;| 4^_p^^d^^ _^,^^^ |x;)) 

nr.,.^r i('e*^'('^+("''/*°)'''«'""~''^/'^^''''(""~''-/'»)(t;|4 o,d . \v) + {v\d . d^ , \v)\ 

is,ss7^io,J- 1 \ I r,-/y,s, rp^-p ,sp. ' ' \ I rp.-p ,sp^ ri-p,Si' 'I 

ra,So P 

+Fi(r„) (e'(x+(-B/*o)(p;(xo-pL/2)-p,(.„-p./2))) ^^1 4^_^^^d^^ _^, ^^ |„)) 

To, So 

where 

>C(x,ro,So)(r,.);(r',.') = Gs^s'{r - p',v' -p,t = Q" )e-' (x+(-S/*o)(P«'(-'-P^72)-p„ (.'-p./2))) (^ _ S^,, ^,,y^^,^,,^)6,^..,) 

-Gs^s' (r - p', r' - p',t = 0+)(l - (5(r',s');(ro,s„)<5s„,i) ■ 



